

JV91 

An Experimental Study of 
Nonlinear Dynamic System Identification 

Greselda I. Stry 1 and D. Joseph Mook 2 

Department of Mechanical and Aerospace Engineering 
State University of New York at Buffalo 
Buffalo, New York 14260 
716-636-3058 

Abstract 

A technique based on the Minimum Model Error optimal estimation approach is 
employed for robust identification of a nonlinear dynamic system. A simple harmonic 
oscillator with quadratic position feedback was simulated on an analog computer. With 
the aid of analog measurements and an assumed linear model, the Minimum Model Error 
Algorithm accurately identifies the quadratic nonlinearity. The tests demonstrate that the 
method is robust with respect to prior ignorance of the nonlinear system model, and with 
respect to measurement record length and regardless of initial conditions. 
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Introduction 

The widespread existence of nonlinear behavior in many dynamic systems is well- 
documented, e.g, Thompson and Stewart [1]; Nayfeh and Mook [2]. In particular, 
virtually every problem associated with orbit estimation, flight trajectory estimation, 
spacecraft dynamics, etc., is known to exhibit nonlinear behavior. Many excellent 
methods for analyzing nonlinear system models have been developed. However, a key 
practical link is often overlooked, namely: How does one obtain an accurate mathematical 
model for the dynamics of a particular complicated nonlinear system? Identification, the 
process of developing an accurate system model from system output measurements, may 
provide the answer. 

Nonlinear systems are commonly described using linear models. Many efficient al- 
gorithms for the identification of linear systems exist and their accuracy and ease of 
application encourages their use. However, linearization does not woric in every appli- 
cation, and even when it does provide a reasonable approximation, the approximation is 
normally limited to a small region about the operating point of linearization. In the case 
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of systems with severe nonlinear behavior using a linear model to describe such a system 
leads to inconsistencies ranging from inaccurate numerical results to misrepresentation 
of the system’s qualitative behavior. Since nonlinearities are seldomly easily charac- 
terized, nonlinear identification techniques may prove beneficial in developing accurate 
mathematical representations of nonlinear systems. 

Numerous methods for the identification of nonlinear systems have been developed 
in the past two decades (Natke, Juang and Gawronski [3]). Among the most widely 
used nonlinear identification methods are describing the nonlinear system using a linear 
model, or representing the nonlinear system in a series expansion and obtaining the 
respective coefficients either by using a regression estimation technique, by minimizing 
a cost functional, by using correlation techniques, or by some other approach. Some 
algorithms use the direct equation approach, while others obtain a graphical representation 
of the nonlinear term(s) and then find an analytical model for the nonlinearity. The 
interested reader can find more information on these nonlinear identification techniques 
in Stry and Mook [4]. 

The diversity of nonlinear identification techniques prompts the choice of an algorithm 
based on the needs of the particular application. Typical criteria to look for are: iterations 
required, robustness in the presence of measurement noise, number of measurements 
needed, robustness with respect to knowledge of the inital conditions, and robustness 
with respect to initial assumptions regarding the form of the nonlinearity. The results 
presented in this paper confirm that the Minimum Model Error algorithm excels in the 
above mentioned requirements. 

In previous papers, the Minimum Model Error algorithm (MME) was explained in 
detail (Mook and Junkins [5]), modified for nonlinear identification (Mook [6]), and 
shown to accurately identify exotic nonlinearities in higher order systems (Stry and 
Mook [4]). In this paper, it is shown how the MME algorithm successfully identifies 
nonlinearities using experimental data. An analytical model representing a harmonic 
oscillator with quadratic position feedback is studied. Output data is obtained from an 
analog computer simulation of the nonlinear system and the quadratic term is accurately 
identified. It is shown that the Minimum Model Error algorithm is capable of identifying 
a nonlinear model which accurately reproduces the analog output regardless of knowledge 
an initially assumed model, initial conditions or record length. 


MME Algorithm 

In this section, we briefly review the MME algorithm and how it is used to identify 
nonlinear dynamic systems. A more detailed explanation may be found in Mook and 
Junkins [5], Mook [6], and Stry and Mook [4], 

The MME may be summarized as follows. Suppose there is a nonlinear system 
whose exact analytical representation is unknown, but for which output measurements 
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are available. Using “normal” means (analysis, finite elements, etc.), a system model, 
denoted assumed model, is constructed. The MME combines the assumed model with 
the measurements to determine the correct form of the nonlinear system. The approach 
consists of adding the to-be-determined correction term to the assumed model. A cost 
functional composed of the weighted sum of the correction term plus measurements 
minus estimate residuals is minimized. The minimization yields optimal state trajectories 
in addition to the correction term. A least squares fit of the correction term is performed 
to find the form of the dynamic model error. 

Consider a forced nonlinear dynamic system which may be modeled in state-space 
form by the equation 

x{t) = A^(t) + £{t) + ( 1 ) 

where x{t) is the n x 1 state vector consisting of the system states, A is the n x n state 
matrix, F(i) is an n x 1 vector of known external excitation, and f[x{t),x{t)) is an 
n x 1 vector which includes all of the system nonlinearities. State-observable discrete 
time domain measurements are available for this system in the form 

y(tfc) = £*(*(**)>**) + £*> *o 

where y{t k ) is an m x 1 measurement vector at time i k , g k is the accurate model of 
the measurement process, and v k represents measurement noise. v k is assumed to be 
a zero-mean, gaussian distributed process of known covariance R k - The measurement 
vector y(t k ) may contain one or more of the system states. To implement MME, assume 
that a model, which is generally not the true system model because of the difficulties 
inherent in obtaining the true system model, is constructed in state-vector form as 

i{t) = Ax(t ) +£(t) (**) 

Here, we show a linear model because in practice, linearization is the most common 
approach to modeling nonlinear systems. MME uses the assumed linear model in (3) 
and the noisy measurements in (2) to find the model error. 

The model error, which might include linear terms as well as unknown nonlinear 
term(s), is represented by the addition of a correction term to the assumed linear model 

as / * \ 

x(t) = Ax{t) + E(t) + i(t) 

where d(t) is the n x 1 correction term (dynamic model error) to be estimated later. 

A cost functional, J, that consists of the weighted integral square of the correction 
term plus the weighted sum square of the measurement-minus-estimated measurement 
residuals, is formed: 



( 5 ) 


+ I* d{r) T Wd{r)dr 

Jto 

where M is the number of measurement times, x(t k ) is the estimated state vector and 
W is a weight matrix to be determined. 

J is minimized with respect to the coirection term, d{t). The necessary conditions 
for the minimization lead to the following two point boundary value problem (TPBVP) 
(see Geering [7]), 


i(t) = M‘) + E(t) + i(t) ( 5 „) 

i(o = -* T m (si,) 

m = -\wm (5c) 

= Uh ) + 2H kRk V y( i k) - 9 k (i(h),tk)] (5 d) 

Hk = fe 

£(*°) = £o or A (to) = 0 

x{t f ) = x f or X(t f ) = 0 ( 5 /) 


where A(<) is a vector of costates (Lagrange multipliers). Estimates of the states and 
of the dynamic model error are produced by the solution of this two-point boundary 
value problem. The estimates depend on the particular value of W. The solution is 
repeated until a value of W is obtained which produces state estimates which satisfy the 
covan ance constraint ”, explained next. 

According to the covariance constraint, the measurement-minus-estimated measure- 
ment residual covariance matrix must match the measurement-minus-truth error covari- 
ance matrix. This may be written as 

(£(<*) - 9_ k m k )^t k )\ T \y(t k ) - g k (x(t k ),t k )} as R k ( 6 ) 

During the minimization, the weight W is varied until the state estimates satisfy the 
covariance constraint, i.e., the left hand side of Eq. (6) is approximately equal to 
the right hand side. The correction term or model error is, therefore, the minimum 
adjustment to the model required for the estimated states to predict the measurements 
with approximately the same covariance as the measurement error. 

After W has been determined such that the state estimates satisfy the covariance 
constraint, the final step in the identification procedure is to use a least squares algorithm 
to fit the model error d(t) to the unknown dynamic teim(s). The error is expanded into 
some combination of linear and nonlinear terms, for example, 

<*(0 = «*(0 + /3* 2 (<) + 7* 3 (f) + . . . (7) 

where a, f3, 7, ... are unknown coefficients to be determined by least squares. The least 
squares approach is explained in detail in Mook[6] 
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The TPBVP represented by Eqs. (5a) to (5f) contains jumps in the costates and, 
consequently, in the correction term. As evident from Eq. (5d), the size of the jump is 
directly proportional to the measurement residual at each measurement time. The noisier 
the measurements, the larger the jump size. A multiple shooting algorithm, developed by 
Mook and Lew [8], converts this jump-discontinuous TPBVP into a set of linear algebraic 
equations which may be solved using any linear equation solver. Multiple shooting also 
facilitates the analysis of a large number of measurements, by processing the solution at 
the end of every set of jumps. 

The multiple shooting algorithm presented by Mook and Lew [8] was used to obtain 
the MME solutions used in the tests presented in this paper. It was assumed in the 
examples that MME obtained the dynamic error term without knowledge of the boundary 
conditions on x, so some distortion of the correction term at the initial and final times 
was expected due to the constraints of Eqs. (5e-5f), i.e., by assuming no state knowledge 
is available at to or t /, we constrain A(t 0 ) = 0 and A (if) = 0. Therefore, in ah test 
cases, the initial and final ten percent of the correction term data was ignored in the 
least squares fit. 


Application Examples 

Two nonlinear equations of motion were studied, which represent the motion of 
an undamped harmonic oscillator with different amounts of quadratic position feedback 
(identical equations may arise in other physical systems as well). The equations in state 
space form are 


(i)-(-. i)(;)*U£uO '*> 

where x is position, and the dot indicates differentiation with respect to time. No forcing 
was applied. 

In the following discussion, Eq. (8) is denoted Model A and Eq. (9) is denoted 
Model B. Different initial conditions were used for each system, for a total of four 
different tests. These are shown in Table 1. 
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Table 1. List of conditions used for each test 


Test # 

x(0) 

v(0) 

x*x 

A1 

0.000 

0.261 

-0.526 

A2 

0.000 

0.523 

-0.526 

B1 

0.000 

0.087 

-1.137 

B2 

0.000 

0.261 

-1.137 


To utilize MME, the linear part of Eqs. (8) and (9) was chosen as the assumed model, 
rendering the model error equivalent to the nonlinear term, c * x 2 . Measurements for the 
MME nonlinear identification were generated by simulating test A1 to B2 on an analog 
computer. Position measurements for all four tests were recorded and nonlinear models 
identified. The results were compared with the analytical position and analytical error 
term data, c * x 2 , which were generated for Models A and B using a digital computer. 
MME proved capable of accurately identifying the nonlinear quadratic term in spite of 
ignorance of the assumed model, true initial conditions and record length. 

Analog computer results 

One hundred position measurements were generated on an EAI-2000 analog com- 
puter for all four test cases. All measurements with a sampling rate of 4 Hertz were 
used in the analysis. Position, velocity, and position squared were chosen as the basis 
functions for the least squares fit. It was uncertain if the analog computer would add 
some damping to the system or if it was able to correctly reproduce the stiffness term. 
By including position and velocity in the least squares fit, stiffness and damping could 
be identified if they existed. The identification procedure yielded the numerical values 
shown in Table 2. 

Table 2. Least Square estimates of the nonlinear terms using measurements 
generated by the analog computer. 


Test 

# 

True 

X 

MME 

X 

True 

V 

MME 

V 

■9 


A1 

0.000 

0.003 

0.000 

0.000 

-0.526 

-0.528 

A2 

0.000 

0.003 

0.000 

0.000 

-0.526 

-0.526 

B1 

0.000 

0.008 

0.000 

0.005 

-1.137 

-1.141 

B2 

0.000 

0.003 

0.000 

0.000 

-1,.137 

-1.135 
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The numerical results for the least squares fit of the error term matched the analytically 
predicted coefficients with great accuracy. Figures (la-4a) show the analytical position, 
analog measurements and position predicted by the MME analysis for all analog tests. 
Figures (lb-4b) show the analytical correction term and the error term estimated by MME. 
In all cases the MME identification produced good state estimates. 

The MME algorithm could accurately identify a nonlinear model regardless of the 
initial conditions. As seen from Figures (la) and (4a) (test A1 and B2), the measured 
position and the analytical position differ significantly. The analytical position was 
digitally recalculated for test A1 and B2 using the initial analog measurements as initial 
conditions instead of the initial conditions presented in Table 1. The results are shown 
in Figures (5a) and (6a). In this set of plots the analytical position and the measurements 
are almost identical. Also, as shown in Figures (5b) and (6b), the analytic correction 
term is much more similar to the estimated correction term, confirming that MME does 
not need any knowledge of the initial or final state vector valued 

MME could identify the nonlinear term accurately idenpendent of the record length. 
In test B1 only 40 measurements were employed in the analysis because subsequent 
measurements were saturated. The nonlinear term is identified very well. 

Note that the data appears to be noiseless, as shown in Figures (la-4a). Successfull 
analysis of noisy data using the MME algorithm can be found in Mook[6] and Stry and 
Mook[4]. 




Figure 1. Test A1 a) Analytical, measured (+), and MME estimated 
O position, b) Analytical and MME estimated f) correction term. 

The MME estimates are essentially identical to the measurements. 

3 It was shown in Eqs. (5e) and (5f), that by setting the initial and final costate values to zero, MME 
does not need any knowledge of the initial or end conditions. 
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Figure 2. Test A2 a) Analytical, measured (+), and MME estimated 
( ) position, b) Analytical and MME estimated (“) correction term. 
The MME estimates are essentially identical to the measurements. 



Figure 3. Test B1 a) Analytical, measured (+), and MME estimated 
( ) position, b) Analytical and MME estimated (“) correction term. 
The MME estimates are essentially identical to the measurements. 
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Figure 4. Test B2 a) Analytical, measured (+), and MME estimated 
(*) position, b) Analytical and MME estimated ( ) correction term. 
The MME estimates are essentially identical to the measurements. 



Figure 5. Test A1 a) Analytical, measured (+), and MME estimated ( ) 
position, b) Analytical and MME estimated (') correction term. The 
analytical position was calculated using as initial conditions the initial 
position and velocity measurements from the analog computer. 
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Figure 6. Test B2 a) Analytical, measured (+), and MME estimated Q 
position, b) Analytical and MME estimated Q correction term. The 
analytical position was calculated using as initial conditions the initial 
position and velocity measurements from the analog computer. 

Conclusion 

In this paper, an MME based algorithm was used to accurately identify the quadratic 
term of a nonlinear harmonic oscillator. Data was obtained from an analog computer 
simulation of the nonlinear system. It is demonstrated that the method is robust with 
respect to (lack of) a priori knowledge of the system dynamics. The identification was 
accurate regardless of initial conditions or data record length. 
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